home *** CD-ROM | disk | FTP | other *** search
/ El Mac 9 / El Mac 9.iso / Shareware / Demos / Igor Demo Pro / 1 PutContentsIn Igor Pro Folder / WaveMetrics Procedures / Image and Contour Plots / Image MagPhase < prev    next >
Encoding:
Text File  |  1996-01-29  |  5.3 KB  |  159 lines  |  [TEXT/IGR0]
  1. // This package performs 2D FFTs on matrix data and shows the results 
  2. // as separate magnitude and phase images or as a combined mag/phase image.
  3. // It adds ImageMagPhase to the Macros menu.
  4. //
  5. // In the dialog from ImageMagPhase:
  6. // Normally, a magnitude only plot is approprate because phase is hard
  7. // to interpret. The combined mag/phase image shows the magnitude as
  8. // intensity and the phase as color (red-green)
  9. // If you don't choose zero in center and your input data is real, the zero frequency
  10. // in the rows dimension will be at the left edge of the image.
  11.  
  12.  
  13. #include <Autosize Images>
  14.  
  15. #pragma rtGlobals= 1
  16.  
  17.  
  18. Macro  ImageMagPhase(w,zeroInCenter,magScaleMode,resultMode)
  19.     String w=  StrVarOrDefault("root:Packages:WMImagesMagPhase:wNameSav","")
  20.     Prompt w,"matrix wave:", Popup WaveList("*",";","")
  21.     Variable zeroInCenter=  NumVarOrDefault("root:Packages:WMImagesMagPhase:zicSav",0)+1
  22.     Prompt zeroInCenter,"zero frequency in center:", Popup "No;Yes"
  23.     Variable magScaleMode= NumVarOrDefault("root:Packages:WMImagesMagPhase:msmSav",2)+1
  24.     Prompt magScaleMode,"Magnitude:", Popup "Linear;Sqrt;Log"
  25.     Variable resultMode= NumVarOrDefault("root:Packages:WMImagesMagPhase:rmSav",2)+1
  26.     Prompt resultMode,"Result Mode:", Popup "Magnitude;Mag and Phase;Combined MagPhase;All"
  27.     
  28.      ImageMagPhaseCombined($w,zeroInCenter-1,magScaleMode-1,resultMode-1)
  29. end
  30.  
  31. //    resultMode:
  32. //    0 -> magnitude
  33. //    1 -> magnitude and phase in different images
  34. //    2 -> combinded mag and phase in one image where mag controls the intensity
  35. //        and phase controls the color (red-green)
  36. //    3 -> all of the above
  37. //    Note: if you choose zeroInCenter, you may note the phase is not symmetrical. This is
  38. //    because the negative frequencies are the complex conj of the positive frequencies
  39. //    If your original data is complex, zero will automatically be in the center
  40. //    Note: in all cases the dc component (zero freq in x and y) is zeroed so it does not
  41. //    overwhelm the rest of the data.
  42. Function ImageMagPhaseCombined(w,zeroInCenter,magScaleMode,resultMode)
  43.     Wave w
  44.     Variable zeroInCenter        // set true if you want zero frequence to be in center (uses more mem)
  45.     Variable magScaleMode        // use zero for linear, 1 for sqrt(), 2 for log
  46.     Variable resultMode            // see above comment
  47.  
  48.  
  49.     // Remember input for next time
  50.     String dfSav= GetDataFolder(1);
  51.     NewDataFolder/O/S root:Packages
  52.     NewDataFolder/O/S WMImagesMagPhase
  53.     String/G wNameSav= NameOfWave(w)
  54.     Variable/G zicSav= zeroInCenter
  55.     Variable/G msmSav= magScaleMode
  56.     Variable/G rmSav= resultMode
  57.     SetDataFolder dfSav
  58.     
  59.     Variable wantMag=  (resultMode==0) %| (resultMode==1) %| (resultMode==3)
  60.     Variable wantPhase= (resultMode==1) %| (resultMode==3)
  61.     Variable wantCombined= (resultMode==2) %| (resultMode==3)
  62.     
  63.     String magName= NameOfWave(w)+"_Mag"
  64.     if( wantMag )
  65.         magName= "::"+magName        // make mag outside tmp data folder
  66.     endif
  67.     String PhaseName= NameOfWave(w)+"_Phase"
  68.     if( wantPhase )
  69.         PhaseName= "::"+PhaseName    // make phase outside tmp data folder
  70.     endif
  71.     NewDataFolder/O/S IMPtmp        // give tmp DF a unique name in case we need to leave it around for debug
  72.  
  73.     Duplicate/O w,$magName
  74.     Wave rmag= $magName
  75.     
  76.     if( (DimSize(w,0) %& 1) )
  77.         zeroInCenter= 1;                // force data to complex if num rows is odd. (requirement of fft)
  78.     endif
  79.  
  80.     if( WaveType(w) %& 4 )            // original data doubles?
  81.         if( zeroInCenter )
  82.             Redimension/C rmag        // C to make complex for zero in cent
  83.         endif
  84.     else
  85.         if( zeroInCenter )
  86.             Redimension/S/C rmag        // S in case it was an integer, C to make complex for zero in cent
  87.         else
  88.             Redimension/S rmag        // S in case it was an integer
  89.         endif
  90.     endif
  91.     FFT rmag                            // rmag may or may not be real
  92.     Wave/C cmag= rmag                // at times mag will be real and complex; pick one
  93.     
  94.     // the following creates separate real valued mag and phase waves
  95.     cmag= r2polar(cmag)
  96.     Duplicate/O cmag,$PhaseName
  97.     Wave rPhase= $PhaseName
  98.     Redimension/R rPhase
  99.     rPhase= imag(cmag)
  100.     Redimension/R cmag
  101.  
  102.     // this zeros the dc component of the mag
  103.     Variable p0= x2pnt(rmag,0)
  104.     Variable q0= (0 - DimOffset(rmag, 1))/DimDelta(rmag,1)    // equiv of y2pnt
  105.     rmag[p0][q0]= 0
  106.  
  107.     if( magScaleMode==1 )
  108.         rmag= sqrt(rmag)
  109.     else
  110.         if( magScaleMode==2 )
  111.             rmag= log(rmag)
  112.         endif
  113.     endif
  114.     
  115.     if( wantCombined )
  116.         String resultWave= NameOfWave(w)+"_MPC"
  117.         String cindexWave= NameOfWave(w)+"_CIndex"
  118.         
  119.         Duplicate/O cmag,::$resultWave
  120.         Wave rw= ::$resultWave
  121.         Redimension/B/U rw
  122.     
  123.         Variable nmags= 20,nphases=8
  124.         Make/N=(20*nphases,3)/W/U/O ::$cindexWave
  125.         Wave mpColors= ::$cindexWave
  126.         mpColors[][0]=65535*(mod(p,nphases)/(nphases-1))*floor(p/nphases)/(nmags-1)
  127.         mpColors[][1]=65535*(((nphases-1)-mod(p,nphases))/(nphases-1))*floor(p/nphases)/(nmags-1)
  128.         mpColors[][2]=0
  129.         
  130.         WaveStats/Q rmag
  131.         rw= round(nmags*((rmag-V_min)/(V_max-V_min)))*nphases +   floor(0.001+(rPhase+pi)*(nphases-1)/(2*pi))
  132.         CheckDisplayed/A rw
  133.         if( V_Flag == 0 )
  134.             Display as "Combined Mag/Phase";AppendImage rw;ModifyImage $resultWave cindex= mpColors
  135.             DoAutoSizeImage(0,1)
  136.             AutoPositionWindow/E
  137.         endif
  138.     endif
  139.     if( wantMag )
  140.         CheckDisplayed/A rmag
  141.         if( V_Flag == 0 )
  142.             Display as "Magnitude";AppendImage rmag;ModifyImage $NameOfWave(rmag)  ctab= {*,*,Grays}
  143.             DoAutoSizeImage(0,1)
  144.             AutoPositionWindow/E
  145.         endif
  146.     endif
  147.     if( wantPhase )
  148.         CheckDisplayed/A rPhase
  149.         if( V_Flag == 0 )
  150.             Display as "Phase";AppendImage rPhase;ModifyImage $NameOfWave(rPhase)  ctab= {*,*,RedWhiteBlue}
  151.             DoAutoSizeImage(0,1)
  152.             AutoPositionWindow/E
  153.         endif
  154.     endif
  155.     KillDataFolder :
  156. end
  157.  
  158.  
  159.